The lattice Boltzmann (LB) method is a powerful approach for hydrodynamics.^{[1, 2]} The essence of the LB method is an intuitively parallel collision-streaming algorithm with discretized position , time t, and microscopic velocity

The ED discretization has been investigated in-depth and a lot of excellent theories have been proposed, such as the small-Mach-number approximation,^[5] the Hermite polynomial expansion,^[6] and the entropic LB model.^[7] According to the Hermite polynomial expansion,^{[6, 8, 9]} for nth-moment-order ED discretization which restores uⁿ moment integral, its can be expressed as

An available smallest quadrature for 2nth QES is the (n+1)th Gauss–Hermite quadrature, which are the zeros of (n+1)th Hermite polynomial. The issue is that the zeros of an Hermite polynomial with degree above 3 cannot fit into nodes, which means that they cannot be expressed as , where v and c stand for integer-valued discrete micro velocity and real-valued lattice constant , respectively. This leads to an off-node lattices in ED discretization. Hence, to construct an on-node lattices for ED discretization—i.e., -type quadrature—one has to manually solve QES, which involves both and . To simplify the notation, in the discussion hereinafter, “lattices” would directly denote “on-node lattices” unless otherwise stated. In practice, a symmetric discrete velocity set is predefined. This avoids the computation of QE with odd exponent k, which significantly simplifies QES, and makes QES purely consist of c and . Employing the skills in Refs. [9] and [10] to deal with QES, a univariate polynomial equation for lattice constant c can be obtained, which separates the co-solving of c and . This leads us to a performable construction of on-node lattices. Actually, this univariate polynomial equation can be directly obtained through a mathematical tool and avoids the tedious QES solving.

In this paper, the partial Gauss–Hermite quadrature (pGHQ) mathematical tool is proposed. The pGHQ is a quadrature rule derived from the Gauss–Hermite quadrature. It keeps the most desirable characteristic of the Gauss–Hermite quadratur; i.e., its quadrature is constructed directly on abscissa polynomial avoiding the co-solving of and in QES. Meanwhile, it offers a performable approach for on-node lattices construction. The on-node lattices construction in the pGHQ scheme is extremely concise. Once a discrete velocity set has been given, a full-range univariate polynomial equation system of its lattice constant c would be directly obtained through pGHQ. Compared with the existing schemes, our approach has the following advantages: (i) the algorithm is extremely concise; (ii) the procedure to construct the univariate polynomial equations is unified for both symmetric and asymmetric lattices; and (iii) the generated univariate polynomial equation system covers full-range quadrature degree of the given . We will elaborate these points in detail in the following.

The theory of pGHQ can be stated as follows: for a q-point abscissa set , whose abscissa polynomial satisfies the orthogonal relationship

Now, we employ pGHQ to construct the univariate polynomial equation of c. The univariate polynomial equation of c essentially is a relation between c and the quadrature degree of the corresponding abscissa set . Once the equation is satisfied, its corresponding possesses a certain quadrature degree, fulfilling QES with a specific order. In classical approaches,^{[8, 10]} it is obtained through manually computing the QES, which needs to construct the QES and separate the co-solving of c and . Now, as the previous discussion shows that HCES is an equation system equivalent to QES but without involving , this relation can be directly constructed by calculating its Hermite polynomial coefficients in abscissa polynomial. Given a predefined with an unknown lattice constant c, we substitute it into the abscissa polynomial with relation , and expand the product

In practice, given a discrete velocity set , the quadrature degree of is required to be as high as possible so that it can be used to construct higher moment degree ED discretization. Therefore, one can start with solving its HCES, where is theoretically the largest. Once this HCES has no real-valued solutions for c, one decreases K by 1. As the construction of HCES shows, this decreasing is actually loosing the constraints on lattice constant c by reducing the governing equations . This procedure is repeated until a real-valued c is found. Its corresponding K gives the quadrature degree of , q + K, which indicates that this set can be used to construct the ED discretization. The construction of is illustrated in Eq. (2). We designate this approach as the pGHQ scheme. It should be noted that there is no limitation on the given discrete velocity set. Given any kind of discrete velocity set, whether it is symmetric or asymmetric, the coefficients can always be obtained and their procedures are unified with the same formulas Eqs. (10)–(12), which is another great advantage of the pGHQ scheme. Hence, the pGHQ scheme supports constructing all kinds of lattices, symmetric or asymmetric.

Compared with the classical approaches,^{[9, 10]} the construction of univariate polynomial equation for lattice constant c in the pGHQ scheme is systematical and general, supporting symmetric and asymmetric lattices and covering all quadrature degrees. The procedure is concise without involving co-solving of c and . It can be mathematically proven that the univariate polynomial equation of c in Refs. [9] and [10] equals HCES. Here, a justification for the Shan scheme^[9] is offered in Section 3 of the Supplementary Material. It is also worth noting that though pGHQ is proposed for the Hermite polynomial expansion theory, it can also be extended into entropic LB model. Actually, it is the mathematical mechanism of a popular entropic LB discretization, the Karlin–Asinari scheme.^[11] The detailed justification can be found in Section 4 of the Supplementary Material. This explains the interesting question^[9]—why, for a given discrete velocity set, does one get the same lattice constant and weights under different schemes, even under different theories?

Since the pGHQ scheme offers a series of HCES covering the full-range quadrature degree and supports all kinds of lattices, a direct application is to construct optimal lattices, which restores the same moment degree with the smallest discrete velocity set. In terms of the pGHQ scheme, given an nth-order moment degree on-node ED discretization, it is to construct a discrete velocity set with the smallest q, whose HCES has real-valued solution for lattice constant c. The theoretically smallest number for q is (n+1), which indicates that its corresponding abscissa polynomial can be expressed as

(I) set up q, start up with the theoretically smallest number q=n+1;

(II) enumerate all the possible q-point discrete velocity sets on [−m, m]

(III) solve HCES for each enumerated discrete velocity set. Identify the set with real-valued c as feasible lattices;

(IV) all the identified feasible sets are local optimal lattices on [−m, m]. If there is no feasible lattice in the enumerated sets, then increase q by 1, repeat steps (II)–(IV).

We find that all of these local optimal lattices keep the symmetric form, . The local optimal abscissa number q on has the relationship with the moment degree n. To verify the feasibility of asymmetric lattices, we continue our search with an extra point. The search shows that for a given n moment degree ED discretization, the available lattices are extremely abundant. Taking n = 3 moment degree ED discretization as an instance, in the range , there are 20 5-point lattices (local optimal lattices) whose discrete velocity set has the form and 34636 6-point lattices, where most of them are asymmetric lattices. Table 1 lists the detailed statistics of our search. As a detailed illustration of the local optimal lattices, Table 2 presents the most compact local optimal lattices on for moment degree ED discretization, whose discrete velocity is as close as possible to 0. To give a specific display of the abundance of the available lattices, Tables 3 and 4 list several symmetric and asymmetric lattices of n = 3 moment degree ED discretization.

Table 1.

Table 1.

Table 1. Statistics of the available on-node lattices for {n = 3,4,5,6,7} moment degree ED discretization on the interval [−10,10]. The local optimal q, the total number of available local optimal lattices, and the total number of available (q+1)-point lattices are respectively listed in columns 2, 3, and 4 of the table. . Moment Local q-point (q+1)-point degree un optimal q lattices lattices 3 5 20 34636 4 7 120 138715 5 9 112 244218 6 11 252 211863 7 13 112 82684

Table 1. Statistics of the available on-node lattices for {n = 3,4,5,6,7} moment degree ED discretization on the interval [−10,10]. The local optimal q, the total number of available local optimal lattices, and the total number of available (q+1)-point lattices are respectively listed in columns 2, 3, and 4 of the table. .

Table 2.

Table 2.

Table 2. The most compact local optimal on-node lattices on [−10,10] for {n = 3,4,5,6,7} moment degree ED discretization. . Moment degree un Lattices Lattice constant c weights 3 {0,±1,±3} 1.1664 {6.3665×10−1, 1.8141×10−1, 2.6196×10−4} 5.5343×10−1 {7.4464×10−2, 4.1859×10−1, 4.4182×10−2} 4 {0,±1,±2,±3} 8.4639×10−1 {4.7667×10−1, 2.3391×10−1, 2.6938×10−2, 8.1213×10−4} 5 {0,±1,±2,±3,±5} 8.1321×10−1 {4.5814×10−1, 2.3734×10−1, 3.2325×10−2, 1.2641×10−3, 8.9773×10−7} 4.7942×10−1 {1.6724×10−1, 3.0315×10−1, 5.3303×10−2, 5.7922×10−2, 2.0013×10−3} 6 {0,±1,±2,±3,±4,±5} 6.8590×10−1 {3.8694×10−1, 2.4178×10−1, 5.8922×10−2, 5.6153×10−3, 2.0652×10−4, 3.2745×10−6} 7 {0,±1,±2,±3,±4,±5,±7} 6.6344×10−1 {3.7428×10−1, 2.4105×10−1, 6.4343×10−2, 7.1316×10−3, 3.2523×10−4, 6.6163×10−6, 3.0509×10−9} 4.3240×10−1 {2.0928×10−1, 2.3312×10−1, 9.4051×10−2, 5.6923×10−2, 7.5008×10−3, 3.7006×10−3, 6.0784×10−5}

Table 2. The most compact local optimal on-node lattices on [−10,10] for {n = 3,4,5,6,7} moment degree ED discretization. .

Table 3.

Table 3.

Table 3. All the available local optimal on-node lattices for n = 3 moment degree ED discretization on the interval [−5,5]. All feasible discrete velocity sets keep the symmetric form, . The and share the same weight . Each feasible lattice has two lattice constants. To save space, the two lattice constants c and their corresponding weights are denoted by the symbols and ± . The rational form is kept for comparison with the existing ED discretizations. . c

Table 3. All the available local optimal on-node lattices for n = 3 moment degree ED discretization on the interval [−5,5]. All feasible discrete velocity sets keep the symmetric form, . The and share the same weight . Each feasible lattice has two lattice constants. To save space, the two lattice constants c and their corresponding weights are denoted by the symbols and ± . The rational form is kept for comparison with the existing ED discretizations. .

Table 4.

Table 4.

Table 4. The available asymmetric on-node lattices for n = 3 moment degree ED discretization. . c w0 w1 w2 w3 w4 w5 {−5, −2, −1, 1, 2, 4} 0.381641 0.019568 0.302751 0.094520 0.505439 0.009237 0.068487 {−4, −3, −1, 1, 2, 4} 0.450877 0.016717 0.054744 0.451349 0.323613 0.127736 0.025841 {−3, −1, 0, 1, 2, 4} 0.521696 0.059199 0.366034 0.198867 0.227904 0.138130 0.009866 {−3, −1, 0, 1, 2, 5} 0.553432 0.076212 0.294489 0.352957 0.040448 0.232265 0.003629

Table 4. The available asymmetric on-node lattices for n = 3 moment degree ED discretization. .

A direct implication of lattices abundance is its impact on the positivity of equilibrium distribution; i.e., the range of macro velocity on which all equilibrium distributions remains positive. Because a negative equilibrium distribution violates the physical nature of particle kinetics, the positivity is considered as a major factor for the numerical stability of the LB method. We have analyzed lattices {0,±1}, {0,±2, ±5}, {0,±1,±3}, {0,±1,±2,±3}, {0,±1,±2,±3,±5}, and {0,±1,±2,±3,±4,±5}. For lattices {0,±2, ±5}, {0,±1,±3}, {0,±1,±2,±3,±5} have two feasible lattice constants c, we take the c with a wider positivity. The analysis shows that lattice {0,±2, ±5} has the widest positivity though its retained moment degree is only u³. Meanwhile, the positivity of the highest retaining-moment-degree lattice {0,±1,±2,±3,±4,±5} is merely better than {0,±1} and {0,±1,±2,±3}. To demonstrate it, Figure 1 plots their equilibrium distributions which firstly go negative as the macro velocity increases. The asymmetric lattice also demonstrates its capability on modifying the positivity on a specific range of U. Figure 2 plots a comparison of lattices {−5, −2, −1, 1, 2, 4} and {0, ±2, ±5}. This shows that the lattice {−5, −2, −1, 1, 2, 4} shifts the positivity range of {0, ±2, ±5} left with approximatively −0.5 on the U-axis. This analysis indicates that the discrete velocity set could be a significant impact to the numerical stability of LB method. It offers a direction to improve LB numerical stability. Our identified lattices also offer a database for further study. However, a detailed investigation is beyond the scope of this paper and shall be addressed in a separate publication.

Fig. 1.

Figure Option ViewDownloadNew Window

Fig. 1. The profiles of first-going-negative equilibrium distributions as a function of U. The figure only renders the positive U-axis. Since all lattices in the figure are symmetric, the positivity on the negative U-axis is the same though the corresponding turns to . The line with symbol is for with and which, as U increases, first goes negative in all equilibrium distributions of is for and c=0.3442 in is for and c=0.5534 in {0,±1,±3}; is for and c=0.8464 in {0,±1,±2,±3}; is for and c=0.4794 in {0,±1,±2,±3,±5}; is for and c=0.6859 in {0,±1,±2,±3,±4,±5}. The inner panel renders their intersections with the U-axis, above which the will become negative. The specific values of intersections for { } are ∼{0.82, 1.70, 1.15, 0.76, 1.25, 0.98}. Since the plotted curves are the first-going-negative equilibrium distributions, then the inner panel demonstrates the lattices positivity range of U.

Fig. 2.

	Figure Option ViewDownloadNew Window
	Fig. 2. The profiles of equilibrium distributions as a function of U: (a) lattice {−5, −2, −1, 1, 2, 4} with a lattice constant c=0.3816; (b) lattice {0, ±2, ±5} with a lattice constant c=0.3442. The label on a curve is its corresponding . The intervals of U between vertical dashed lines are lattices positivity ranges.

We propose a new mathematical tool, pGHQ, to construct on-node LB lattices under a Cartesian coordinate system in this paper. To the best of our knowledge, this is the first time to derive and employ this mathematical tool in the context of the LB method. The pGHQ is general. It can be extended into the entropic LB model, even though it was first proposed for the Hermite polynomial expansion theory. The pGHQ scheme avoids the tedious QES solving. Compared with the existing classical approaches, our scheme has the following advantages: (i) the algorithm is extremely concise; (ii) the procedure of constructing univariate polynomial equations is unified for both symmetric and asymmetric lattices; and (iii) the generated univariate polynomial equation system covers full-range quadrature degree of the given . We employ the pGHQ scheme to search the local optimal and asymmetric lattices on [−10,10] for {n = 3,4,5,6,7} moment degree ED discretization. The search reveals a surprising abundance of available lattices. Our brief analysis shows that the discrete velocity set is significant to the positivity of equilibrium distribution, which is one major impact to the numerical stability of LB method. Hence, the results of the pGHQ scheme lay a foundation for further investigation to improve the numerical stability of the LB method by modifying the discrete velocity set.